home *** CD-ROM | disk | FTP | other *** search
/ Amiga Tools 3 / Amiga Tools 3.iso / grafik / raytracing / rayshade-4.0.6.3 / urt / lib / colorquant.c < prev    next >
Encoding:
C/C++ Source or Header  |  1991-08-09  |  20.2 KB  |  744 lines
  1. /*
  2.  * This software is copyrighted as noted below.  It may be freely copied,
  3.  * modified, and redistributed, provided that the copyright notice is 
  4.  * preserved on all copies.
  5.  * 
  6.  * There is no warranty or other guarantee of fitness for this software,
  7.  * it is provided solely "as is".  Bug reports or fixes may be sent
  8.  * to the author, who may or may not act on them as he desires.
  9.  *
  10.  * You may not include this software in a program or other software product
  11.  * without supplying the source, or without informing the end-user that the 
  12.  * source is available for no extra charge.
  13.  *
  14.  * If you modify this software, you should include a notice giving the
  15.  * name of the person performing the modification, the date of modification,
  16.  * and the reason for such modification.
  17.  */
  18. /*
  19.  * colorquant.c
  20.  *
  21.  * Perform variance-based color quantization on a "full color" image.
  22.  * Author:    Craig Kolb
  23.  *        Department of Mathematics
  24.  *        Yale University
  25.  *        kolb@yale.edu
  26.  * Date:    Tue Aug 22 1989
  27.  * Copyright (C) 1989 Craig E. Kolb
  28.  * $Id: colorquant.c,v 3.0.1.2 90/11/29 15:18:04 spencer Exp $
  29.  *
  30.  * $Log:    colorquant.c,v $
  31.  * Revision 3.0.1.2  90/11/29  15:18:04  spencer
  32.  * Remove a typo.
  33.  * 
  34.  * 
  35.  * Revision 3.0.1.1  90/11/19  16:59:48  spencer
  36.  * Use inv_cmap instead of find_colors.
  37.  * Changes to process multiple files into one colormap (accum_hist arg).
  38.  * Delete 'otherimages' argument -- unnecessary with faster inv_cmap code.
  39.  * 
  40.  * 
  41.  * Revision 3.0  90/08/03  15:20:11  spencer
  42.  * Establish version 3.0 base.
  43.  * 
  44.  * Revision 1.6  90/07/29  08:06:06  spencer
  45.  * If HUGE isn't defined, make it HUGE_VAL (for ansi).
  46.  * 
  47.  * Revision 1.5  90/07/26  17:25:48  rgb
  48.  * Added a parameter to colorquant for rgbmap construction.
  49.  * 
  50.  * Revision 1.4  90/07/13  14:53:31  spencer
  51.  * Get rid of ARB_ARG sh*t.
  52.  * Change a couple of vars to double so that HUGE won't cause problems.
  53.  * 
  54.  * Revision 1.3  90/06/28  21:42:56  spencer
  55.  * Declare internal functions properly.
  56.  * Delete unused global variable.
  57.  * 
  58.  * Revision 1.2  90/06/28  13:18:56  spencer
  59.  * Make internal functions static.
  60.  * Build entire RGB cube, not just those colors used by this image.  This
  61.  * is so dithering will work.
  62.  * 
  63.  * Revision 1.1  90/06/18  20:45:17  spencer
  64.  * Initial revision
  65.  * 
  66.  * Revision 1.3  89/12/03  18:27:16  craig
  67.  * Removed bogus integer casts in distance calculation in makenearest().
  68.  * 
  69.  * Revision 1.2  89/12/03  18:13:12  craig
  70.  * FindCutpoint now returns FALSE if the given box cannot be cut.  This
  71.  * to avoid overflow problems in CutBox.
  72.  * "whichbox" in GreatestVariance() is now initialized to 0.
  73.  * 
  74.  */
  75. #include <math.h>
  76. #include <stdio.h>
  77.  
  78. /* Ansi uses HUGE_VAL. */
  79. #ifndef HUGE
  80. #define HUGE HUGE_VAL
  81. #endif
  82.  
  83. static void QuantHistogram();
  84. static void BoxStats();
  85. static void UpdateFrequencies();
  86. static void ComputeRGBMap();
  87. static void SetRGBmap();
  88. #ifdef slow_color_map
  89. static void find_colors();
  90. static int  getneighbors();
  91. static int  makenearest();
  92. #else
  93. extern void inv_cmap();
  94. #endif
  95. static int  CutBoxes();
  96. static int  CutBox();
  97. static int  GreatestVariance();
  98. static int  FindCutpoint();
  99.  
  100.  
  101.  
  102. /* 
  103.  * The colorquant() routine has been tested on an Iris 4D70 workstation,
  104.  * a Sun 3/60 workstation, and (to some extent), a Macintosh.
  105.  * 
  106.  * Calls to bzero() may have to be replaced with the appropriate thing on
  107.  * your system.  (bzero(ptr, len) writes 'len' 0-bytes starting at the location
  108.  * pointed to by ptr.)
  109.  * 
  110.  * Although I've tried to avoid integer overflow problems where ever possible,
  111.  * it's likely I've missed a spot where an 'int' should really be a 'long'.
  112.  * (On the machine this was developed on, an int == long == 32 bits.)
  113.  * 
  114.  * Note that it's quite easy to optimize this code for a given value for
  115.  * 'bits'.  In addition, if you assume bits is small and
  116.  * that the total number of pixels is relatively small, there are several
  117.  * places that integer arithmetic may be substituted for floating-point.
  118.  * One such place is the loop in BoxStats -- mean and var need not necessary
  119.  * be floats.
  120.  * 
  121.  * As things stand, the maximum number of colors to which an image may
  122.  * be quantized is 256.  This limit may be overcome by changing rgbmap and
  123.  * colormap from pointers to characters to pointers to something larger.
  124.  */
  125.  
  126. /*
  127.  * Maximum number of colormap entries.  To make larger than 2^8, the rgbmap
  128.  * type will have to be changed from unsigned chars to something larger.
  129.  */
  130. #define MAXCOLORS        256
  131. /*
  132.  * Value corrersponding to full intensity in colormap.  The values placed
  133.  * in the colormap are scaled to be between zero and this number.  Note
  134.  * that anything larger than 255 is going to lead to problems, as the
  135.  * colormap is declared as an unsigned char.
  136.  */
  137. #define FULLINTENSITY        255
  138. #define MAX(x,y)    ((x) > (y) ? (x) : (y))
  139.  
  140. /*
  141.  * Readability constants.
  142.  */
  143. #define REDI        0    
  144. #define GREENI        1
  145. #define BLUEI        2    
  146. #define TRUE        1
  147. #define FALSE        0
  148.  
  149. typedef struct {
  150.     float        weightedvar,        /* weighted variance */
  151.             mean[3];        /* centroid */
  152.     unsigned long     weight,            /* # of pixels in box */
  153.             freq[3][MAXCOLORS];    /* Projected frequencies */
  154.     int         low[3], high[3];    /* Box extent */
  155. } Box;
  156.  
  157. static unsigned long    *Histogram,        /* image histogram */
  158.             NPixels,        /* # of pixels in an image*/
  159.             SumPixels;        /* total # of pixels */
  160. static unsigned int    Bits,            /* # significant input bits */
  161.             ColormaxI;        /* # of colors, 2^Bits */
  162. static Box        *Boxes;            /* Array of color boxes. */
  163.  
  164. /*
  165.  * Perform variance-based color quantization on a 24-bit image.
  166.  *
  167.  * Input consists of:
  168.  *    red, green, blue    Arrays of red, green and blue pixel
  169.  *                intensities stored as unsigned characters.
  170.  *                The color of the ith pixel is given
  171.  *                by red[i] green[i] and blue[i].  0 indicates
  172.  *                zero intensity, 255 full intensity.
  173.  *    pixels            The length of the red, green and blue arrays
  174.  *                in bytes, stored as an unsigned long int.
  175.  *    colormap        Points to the colormap.  The colormap
  176.  *                consists of red, green and blue arrays.
  177.  *                The red/green/blue values of the ith
  178.  *                colormap entry are given respectively by
  179.  *                colormap[0][i], colormap[1][i] and
  180.  *                colormap[2][i].  Each entry is an unsigned char.
  181.  *    colors            The number of colormap entries, stored
  182.  *                as an integer.    
  183.  *    bits            The number of significant bits in each entry
  184.  *                of the red, greena and blue arrays. An integer.
  185.  *    rgbmap            An array of unsigned chars of size (2^bits)^3.
  186.  *                This array is used to map from pixels to
  187.  *                colormap entries.  The 'prequantized' red,
  188.  *                green and blue components of a pixel
  189.  *                are used as an index into rgbmap to retrieve
  190.  *                the index which should be used into the colormap
  191.  *                to represent the pixel.  In short:
  192.  *                index = rgbmap[(((r << bits) | g) << bits) | b];
  193.  *     fast            If non-zero, the rgbmap will be constructed
  194.  *                quickly.  If zero, the rgbmap will be built
  195.  *                much slower, but more accurately.  In most
  196.  *                cases, fast should be non-zero, as the error
  197.  *                introduced by the approximation is usually
  198.  *                small.  'Fast' is stored as an integer.
  199.  *    accum_hist        If non-zero the histogram will accumulate and 
  200.  *                reflect pixels from multiple images.
  201.  *                when 1, the histogram will be initialized and
  202.  *                summed, but not thrown away OR processed. when 
  203.  *                2 the image RGB will be added to it.  When 3 
  204.  *                Boxes are cut and a colormap and rgbmap
  205.  *                are be returned, Histogram is freed too.
  206.  *                When zero, all code is executed as per normal.
  207.  *
  208.  * colorquant returns the number of colors to which the image was
  209.  * quantized.
  210.  */
  211. #define INIT_HIST 1
  212. #define USE_HIST 2
  213. #define PROCESS_HIST 3
  214. int
  215. colorquant(red,green,blue,pixels,colormap,colors,bits,rgbmap,fast,accum_hist)
  216. unsigned char *red, *green, *blue;
  217. unsigned long pixels;
  218. unsigned char *colormap[3];
  219. int colors, bits;
  220. unsigned char *rgbmap;
  221. int fast, accum_hist;
  222. {
  223.     int    i,            /* Counter */
  224.     OutColors,            /* # of entries computed */
  225.     Colormax;            /* quantized full-intensity */ 
  226.     float    Cfactor;    /* Conversion factor */
  227.  
  228.     if (accum_hist < 0 || accum_hist > PROCESS_HIST)
  229.     fprintf(stderr, "colorquant: bad value for accum_hist\n");
  230.  
  231.     ColormaxI = 1 << bits;    /* 2 ^ Bits */
  232.     Colormax = ColormaxI - 1;
  233.     Bits = bits;
  234.     NPixels = pixels;
  235.     Cfactor = (float)FULLINTENSITY / Colormax;
  236.  
  237.     if (! accum_hist || accum_hist == INIT_HIST) {
  238.     Histogram = (unsigned long *)calloc(ColormaxI*ColormaxI*ColormaxI, 
  239.                         sizeof(long));
  240.     Boxes = (Box *)malloc(colors * sizeof(Box));
  241.     /*
  242.      * Zero-out the projected frequency arrays of the largest box.
  243.      */
  244.     bzero(Boxes->freq[0], ColormaxI * sizeof(unsigned long));
  245.     bzero(Boxes->freq[1], ColormaxI * sizeof(unsigned long));
  246.     bzero(Boxes->freq[2], ColormaxI * sizeof(unsigned long));
  247.     SumPixels = 0;
  248.     }
  249.  
  250.     SumPixels += NPixels;
  251.  
  252.     if ( accum_hist != PROCESS_HIST ) 
  253.     QuantHistogram(red, green, blue, &Boxes[0]);
  254.  
  255.     if ( !accum_hist || accum_hist == PROCESS_HIST) {
  256.     OutColors = CutBoxes(Boxes, colors);
  257.  
  258.     /*
  259.      * We now know the set of representative colors.  We now
  260.      * must fill in the colormap and convert the representatives
  261.      * from their 'prequantized' range to 0-FULLINTENSITY.
  262.      */
  263.     for (i = 0; i < OutColors; i++) {
  264.         colormap[0][i] =
  265.         (unsigned char)(Boxes[i].mean[REDI] * Cfactor + 0.5);
  266.         colormap[1][i] =
  267.         (unsigned char)(Boxes[i].mean[GREENI] * Cfactor + 0.5);
  268.         colormap[2][i] =
  269.         (unsigned char)(Boxes[i].mean[BLUEI] * Cfactor + 0.5);
  270.     }
  271.  
  272.     ComputeRGBMap(Boxes, OutColors, rgbmap, bits, colormap, fast);
  273.  
  274.     free((char *)Histogram);
  275.     free((char *)Boxes);
  276.     return OutColors;    /* Return # of colormap entries */
  277.     }
  278.     return 0;
  279. }
  280.  
  281. /*
  282.  * Compute the histogram of the image as well as the projected frequency
  283.  * arrays for the first world-encompassing box.
  284.  */
  285. static void
  286. QuantHistogram(r, g, b, box)
  287. register unsigned char *r, *g, *b;
  288. Box *box;
  289. {
  290.     unsigned long *rf, *gf, *bf, i;
  291.  
  292.     rf = box->freq[0];
  293.     gf = box->freq[1];
  294.     bf = box->freq[2];
  295.  
  296.     /*
  297.      * We compute both the histogram and the proj. frequencies of
  298.      * the first box at the same time to save a pass through the
  299.      * entire image. 
  300.      */
  301.  
  302.     for (i = 0; i < NPixels; i++) {
  303.         rf[*r]++;
  304.         gf[*g]++;
  305.         bf[*b]++;
  306.         Histogram[((((*r++)<<Bits)|(*g++))<<Bits)|(*b++)]++;
  307.     }
  308. }
  309.  
  310. /*
  311.  * Interatively cut the boxes.
  312.  */
  313. static
  314. CutBoxes(boxes, colors) 
  315. Box    *boxes;
  316. int    colors;
  317. {
  318.     int curbox;
  319.  
  320.     boxes[0].low[REDI] = boxes[0].low[GREENI] = boxes[0].low[BLUEI] = 0;
  321.     boxes[0].high[REDI] = boxes[0].high[GREENI] =
  322.                   boxes[0].high[BLUEI] = ColormaxI;
  323.     boxes[0].weight = SumPixels;
  324.  
  325.     BoxStats(&boxes[0]);
  326.  
  327.     for (curbox = 1; curbox < colors; curbox++) {
  328.         if (CutBox(&boxes[GreatestVariance(boxes, curbox)],
  329.                &boxes[curbox]) == FALSE)
  330.                 break;
  331.     }
  332.  
  333.     return curbox;
  334. }
  335.  
  336. /*
  337.  * Return the number of the box in 'boxes' with the greatest variance.
  338.  * Restrict the search to those boxes with indices between 0 and n-1.
  339.  */
  340. static
  341. GreatestVariance(boxes, n)
  342. Box *boxes;
  343. int n;
  344. {
  345.     register int i, whichbox = 0;
  346.     float max;
  347.  
  348.     max = -1;
  349.     for (i = 0; i < n; i++) {
  350.         if (boxes[i].weightedvar > max) {
  351.             max = boxes[i].weightedvar;
  352.             whichbox = i;
  353.         }
  354.     }
  355.     return whichbox;
  356. }
  357.  
  358. /*
  359.  * Compute mean and weighted variance of the given box.
  360.  */
  361. static void
  362. BoxStats(box)
  363. register Box *box;
  364. {
  365.     register int i, color;
  366.     unsigned long *freq;
  367.     float mean, var;
  368.  
  369.     if(box->weight == 0) {
  370.         box->weightedvar = 0;
  371.         return;
  372.     }
  373.  
  374.     box->weightedvar = 0.;
  375.     for (color = 0; color < 3; color++) {
  376.         var = mean = 0;
  377.         i = box->low[color];
  378.         freq = &box->freq[color][i];
  379.         for (; i < box->high[color]; i++, freq++) {
  380.             mean += i * *freq;
  381.             var += i*i* *freq;
  382.         }
  383.         box->mean[color] = mean / (float)box->weight;
  384.         box->weightedvar += var - box->mean[color]*box->mean[color]*
  385.                     (float)box->weight;
  386.     }
  387.     box->weightedvar /= SumPixels;
  388. }
  389.  
  390. /*
  391.  * Cut the given box.  Returns TRUE if the box could be cut, FALSE otherwise.
  392.  */
  393. static
  394. CutBox(box, newbox)
  395. Box *box, *newbox;
  396. {
  397.     int i;
  398.     double totalvar[3];
  399.     Box newboxes[3][2];
  400.  
  401.     if (box->weightedvar == 0. || box->weight == 0)
  402.         /*
  403.          * Can't cut this box.
  404.          */
  405.         return FALSE;
  406.  
  407.     /*
  408.      * Find 'optimal' cutpoint along each of the red, green and blue
  409.      * axes.  Sum the variances of the two boxes which would result
  410.      * by making each cut and store the resultant boxes for 
  411.      * (possible) later use.
  412.      */
  413.     for (i = 0; i < 3; i++) {
  414.         if (FindCutpoint(box, i, &newboxes[i][0], &newboxes[i][1]))
  415.             totalvar[i] = newboxes[i][0].weightedvar +
  416.                 newboxes[i][1].weightedvar;
  417.         else
  418.             totalvar[i] = HUGE;
  419.     }
  420.  
  421.     /*
  422.      * Find which of the three cuts minimized the total variance
  423.      * and make that the 'real' cut.
  424.      */
  425.     if (totalvar[REDI] <= totalvar[GREENI] &&
  426.         totalvar[REDI] <= totalvar[BLUEI]) {
  427.         *box = newboxes[REDI][0];
  428.         *newbox = newboxes[REDI][1];
  429.     } else if (totalvar[GREENI] <= totalvar[REDI] &&
  430.          totalvar[GREENI] <= totalvar[BLUEI]) {
  431.         *box = newboxes[GREENI][0];
  432.         *newbox = newboxes[GREENI][1];
  433.     } else {
  434.         *box = newboxes[BLUEI][0];
  435.         *newbox = newboxes[BLUEI][1];
  436.     }
  437.  
  438.     return TRUE;
  439. }
  440.  
  441. /*
  442.  * Compute the 'optimal' cutpoint for the given box along the axis
  443.  * indcated by 'color'.  Store the boxes which result from the cut
  444.  * in newbox1 and newbox2.
  445.  */
  446. static
  447. FindCutpoint(box, color, newbox1, newbox2)
  448. Box *box, *newbox1, *newbox2;
  449. int color;
  450. {
  451.     float u, v, max;
  452.     int i, maxindex, minindex, cutpoint;
  453.     unsigned long optweight, curweight;
  454.  
  455.     if (box->low[color] + 1 == box->high[color])
  456.         return FALSE;    /* Cannot be cut. */
  457.     minindex = (int)((box->low[color] + box->mean[color]) * 0.5);
  458.     maxindex = (int)((box->mean[color] + box->high[color]) * 0.5);
  459.  
  460.     cutpoint = minindex;
  461.     optweight = box->weight;
  462.  
  463.     curweight = 0;
  464.     for (i = box->low[color] ; i < minindex ; i++)
  465.         curweight += box->freq[color][i];
  466.     u = 0.;
  467.     max = -1;
  468.     for (i = minindex; i <= maxindex ; i++) {
  469.         curweight += box->freq[color][i];
  470.         if (curweight == box->weight)
  471.             break;
  472.         u += (float)(i * box->freq[color][i]) /
  473.                     (float)box->weight;
  474.         v = ((float)curweight / (float)(box->weight-curweight)) *
  475.                 (box->mean[color]-u)*(box->mean[color]-u);
  476.         if (v > max) {
  477.             max = v;
  478.             cutpoint = i;
  479.             optweight = curweight;
  480.         }
  481.     }
  482.     cutpoint++;
  483.     *newbox1 = *newbox2 = *box;
  484.     newbox1->weight = optweight;
  485.     newbox2->weight -= optweight;
  486.     newbox1->high[color] = cutpoint;
  487.     newbox2->low[color] = cutpoint;
  488.     UpdateFrequencies(newbox1, newbox2);
  489.     BoxStats(newbox1);
  490.     BoxStats(newbox2);
  491.  
  492.     return TRUE;    /* Found cutpoint. */
  493. }
  494.  
  495. /*
  496.  * Update projected frequency arrays for two boxes which used to be
  497.  * a single box.
  498.  */
  499. static void
  500. UpdateFrequencies(box1, box2)
  501. register Box *box1, *box2;
  502. {
  503.     register unsigned long myfreq, *h;
  504.     register int b, g, r;
  505.     int roff;
  506.  
  507.     bzero(box1->freq[0], ColormaxI * sizeof(unsigned long));
  508.     bzero(box1->freq[1], ColormaxI * sizeof(unsigned long));
  509.     bzero(box1->freq[2], ColormaxI * sizeof(unsigned long)); 
  510.  
  511.     for (r = box1->low[0]; r < box1->high[0]; r++) {
  512.         roff = r << Bits;
  513.         for (g = box1->low[1];g < box1->high[1]; g++) {
  514.             b = box1->low[2];
  515.             h = Histogram + (((roff | g) << Bits) | b);
  516.             for (; b < box1->high[2]; b++) {
  517.                 if ((myfreq = *h++) == 0)
  518.                     continue;
  519.                 box1->freq[0][r] += myfreq;
  520.                 box1->freq[1][g] += myfreq;
  521.                 box1->freq[2][b] += myfreq;
  522.                 box2->freq[0][r] -= myfreq;
  523.                 box2->freq[1][g] -= myfreq;
  524.                 box2->freq[2][b] -= myfreq;
  525.             }
  526.         }
  527.     }
  528. }
  529.  
  530. /*
  531.  * Compute RGB to colormap index map.
  532.  */
  533. static void
  534. ComputeRGBMap(boxes, colors, rgbmap, bits, colormap, fast)
  535. Box *boxes;
  536. int colors;
  537. unsigned char *rgbmap, *colormap[3];
  538. int bits, fast;
  539. {
  540.     register int i;
  541.  
  542.     if (fast) {
  543.         /*
  544.          * The centroid of each box serves as the representative
  545.          * for each color in the box.
  546.          */
  547.         for (i = 0; i < colors; i++)
  548.             SetRGBmap(i, &boxes[i], rgbmap, bits);
  549.     } else
  550.         /*
  551.          * Find the 'nearest' representative for each
  552.          * pixel.
  553.          */
  554. #ifdef slow_color_map
  555.         find_colors(boxes, colors, rgbmap, bits, colormap, 0);
  556. #else
  557.         inv_cmap(colors, colormap, bits, Histogram, rgbmap);
  558. #endif
  559. }
  560.  
  561. /*
  562.  * Make the centroid of "boxnum" serve as the representative for
  563.  * each color in the box.
  564.  */
  565. static void
  566. SetRGBmap(boxnum, box, rgbmap, bits)
  567. int boxnum;
  568. Box *box;
  569. unsigned char *rgbmap;
  570. int bits;
  571. {
  572.     register int r, g, b;
  573.     
  574.     for (r = box->low[REDI]; r < box->high[REDI]; r++) {
  575.         for (g = box->low[GREENI]; g < box->high[GREENI]; g++) {
  576.             for (b = box->low[BLUEI]; b < box->high[BLUEI]; b++) {
  577.                 rgbmap[(((r<<bits)|g)<<bits)|b]=(char)boxnum;
  578.             }
  579.         }
  580.     }
  581. }
  582.  
  583. #ifdef slow_color_map
  584. /*
  585.  * Form colormap and NearestColor arrays.
  586.  */
  587. static void
  588. find_colors(boxes, colors, rgbmap, bits, colormap, otherimages)
  589. int colors;
  590. Box *boxes;
  591. unsigned char *rgbmap;
  592. int bits;
  593. unsigned char *colormap[3];
  594. int otherimages;
  595. {
  596.     register int i;
  597.     int num, *neighbors;
  598.  
  599.     neighbors = (int *)malloc(colors * sizeof(int));
  600.  
  601.     /*
  602.      * Form map of representative (nearest) colors.
  603.      */
  604.     for (i = 0; i < colors; i++) {
  605.         /*
  606.          * Create list of candidate neighbors and
  607.          * find closest representative for each
  608.          * color in the box.
  609.          */
  610.         num = getneighbors(boxes, i, neighbors, colors, colormap);
  611.         makenearest(boxes, i, num, neighbors, rgbmap, bits, colormap, 
  612.                 otherimages);
  613.     }
  614.     free((char *)neighbors);
  615. }
  616.  
  617. /*
  618.  * In order to minimize our search for 'best representative', we form the
  619.  * 'neighbors' array.  This array contains the number of the boxes whose
  620.  * centroids *might* be used as a representative for some color in the
  621.  * current box.  We need only consider those boxes whose centroids are closer
  622.  * to one or more of the current box's corners than is the centroid of the
  623.  * current box. 'Closeness' is measured by Euclidean distance.
  624.  */
  625. static
  626. getneighbors(boxes, num, neighbors, colors, colormap)
  627. Box *boxes;
  628. int num, colors, *neighbors;
  629. unsigned char *colormap[3];
  630. {
  631.     register int i, j;
  632.     Box *bp;
  633.     float dist, LowR, LowG, LowB, HighR, HighG, HighB, ldiff, hdiff;
  634.  
  635.     bp = &boxes[num];
  636.  
  637.     ldiff = bp->low[REDI] - bp->mean[REDI];
  638.     ldiff *= ldiff;
  639.     hdiff = bp->high[REDI] - bp->mean[REDI];
  640.     hdiff *= hdiff;
  641.     dist = MAX(ldiff, hdiff);
  642.  
  643.     ldiff = bp->low[GREENI] - bp->mean[GREENI];
  644.     ldiff *= ldiff;
  645.     hdiff = bp->high[GREENI] - bp->mean[GREENI];
  646.     hdiff *= hdiff;
  647.     dist += MAX(ldiff, hdiff);
  648.  
  649.     ldiff = bp->low[BLUEI] - bp->mean[BLUEI];
  650.     ldiff *= ldiff;
  651.     hdiff = bp->high[BLUEI] - bp->mean[BLUEI];
  652.     hdiff *= hdiff;
  653.     dist += MAX(ldiff, hdiff);
  654.  
  655. #ifdef IRIS
  656.     dist = fsqrt(dist);
  657. #else
  658.     dist = (float)sqrt((double)dist);
  659. #endif
  660.  
  661.     /*
  662.      * Loop over all colors in the colormap, the ith entry of which
  663.      * corresponds to the ith box.
  664.      *
  665.      * If the centroid of a box is as close to any corner of the
  666.      * current box as is the centroid of the current box, add that
  667.      * box to the list of "neighbors" of the current box.
  668.      */
  669.     HighR = (float)bp->high[REDI] + dist;
  670.     HighG = (float)bp->high[GREENI] + dist;
  671.     HighB = (float)bp->high[BLUEI] + dist;
  672.     LowR = (float)bp->low[REDI] - dist;
  673.     LowG = (float)bp->low[GREENI] - dist;
  674.     LowB = (float)bp->low[BLUEI] - dist;
  675.     for (i = j = 0, bp = boxes; i < colors; i++, bp++) {
  676.         if (LowR <= bp->mean[REDI] && HighR >= bp->mean[REDI] &&
  677.             LowG <= bp->mean[GREENI] && HighG >= bp->mean[GREENI] &&
  678.             LowB <= bp->mean[BLUEI] && HighB >= bp->mean[BLUEI])
  679.             neighbors[j++] = i;
  680.     }
  681.  
  682.     return j;    /* Return the number of neighbors found. */
  683. }
  684.  
  685. /*
  686.  * Assign representative colors to every pixel in a given box through
  687.  * the construction of the NearestColor array.  For each color in the
  688.  * given box, we look at the list of neighbors passed to find the
  689.  * one whose centroid is closest to the current color.
  690.  */
  691. static
  692. makenearest(boxes, boxnum, nneighbors, neighbors, rgbmap, bits, colormap, 
  693.         otherimages)
  694. Box *boxes;
  695. int boxnum;
  696. int nneighbors, *neighbors, bits;
  697. unsigned char *rgbmap, *colormap[3];
  698. int otherimages;
  699. {
  700.     register int n, b, g, r;
  701.     double rdist, gdist, bdist, dist, mindist;
  702.     int which, *np;
  703.     Box *box;
  704.  
  705.     box = &boxes[boxnum];
  706.  
  707.     for (r = box->low[REDI]; r < box->high[REDI]; r++) {
  708.         for (g = box->low[GREENI]; g < box->high[GREENI]; g++) {
  709.             for (b = box->low[BLUEI]; b < box->high[BLUEI]; b++) {
  710. /*
  711.  * The following "if" is to avoid doing extra work when the RGBmap is
  712.  * not going to be used for other images.
  713.  */
  714.                     if ((!otherimages) && 
  715.                     (Histogram[(((r<<bits)|g)<<bits)|b] == 0))
  716.                     continue;
  717.  
  718.                 mindist = HUGE;
  719.                 /*
  720.                  * Find the colormap entry which is
  721.                  * closest to the current color.
  722.                  */
  723.                 np = neighbors;
  724.                 for (n = 0; n < nneighbors; n++, np++) {
  725.                     rdist = r-boxes[*np].mean[REDI];
  726.                     gdist = g-boxes[*np].mean[GREENI];
  727.                     bdist = b-boxes[*np].mean[BLUEI];
  728.                     dist = rdist*rdist + gdist*gdist + bdist*bdist;
  729.                     if (dist < mindist) {
  730.                         mindist = dist;
  731.                         which = *np; 
  732.                     }
  733.                 }
  734.                 /*
  735.                  * The colormap entry closest to this
  736.                  * color is used as a representative.
  737.                  */
  738.                 rgbmap[(((r<<bits)|g)<<bits)|b] = which;
  739.             }
  740.         }
  741.     }
  742. }
  743. #endif /* slow_color_map */
  744.